# 加载相关包
library(openxlsx)
library(geepack)
library(mice)
library(reshape2)
# 加载自定义函数
source("../func/gee_run.r")

# 1. 读取数据，数据预处理
data_test = read.xlsx("../input/data918.xlsx", 3)
'
10-09：已知存在问题
# 1.1 删除 6 月的观测、且 8 月只有一个观测的个体
data_tmp = data_test[data_test$month != 6 & data_test$month != 8, ]
extract_id = unique(data_tmp$id)
data_tmp = data_test[data_test$id %in% extract_id, ]
data_tmp = data_tmp[data_tmp$month != 6, ]
'

'10-14：修改逻辑'
# 1.1 保留 8~24 月至少两条观测的个体
data_cast = dcast(data_test, id ~ month)        # 通过 View(data_cast) 查看重铸后的数据
data_cast = transform(data_cast, count = `8` + `12` + `18` + `24`)   # count 计算 8~24 月填了多少次
extract_id = data_cast[data_cast$count > 1, ]$id                     # 包含 8~24 月至少填了两次的个体的向量
data_tmp = data_test[(data_test$id %in% extract_id) & (data_test$month != 6), ]   # 取得这些个体的所有观测（但不要 6 月份的观测）
# 1.2 提取数据，确定自变量、因变量、协变量
data_anal = data_tmp[, c(1:2, 9:13, 14:25, 27, 29:39)]
varsx = colnames(data_test)[c(14:25, 27)]  # 自变量
varsy = colnames(data_test)[c(9:13)]       # 因变量
varsc = colnames(data_test)[c(29:39)]      # 协变量
# 1.3 确保变量类型一致：因变量为数值，自变量、协变量为因子
data_anal$id = as.character(data_anal$id)
data_anal[varsy] = lapply(data_anal[varsy], as.numeric)
data_anal[varsx] = lapply(data_anal[varsx], as.factor)
data_anal[varsc] = lapply(data_anal[varsc], as.factor)


# 2. 填补缺失值并跑 GEE
# 2.1 mice
data_tmp = unique(data_anal[, c("id", varsc)])
miss = md.pattern(data_tmp)                # 查看缺失数据的模式
imp = mice(data_tmp, m = 5, seed = 1)      # 填补
mice_data = complete(imp, action = 1)      # 懒惰的人选择第一个
mice_data = merge(data_anal[, c("id", "month", varsy, varsx)], mice_data)   # 合并
write.xlsx(mice_data, '../output/1009_mice_data.xlsx', overwrite = T)

## 2.2 预处理
mice_data = read.xlsx('../output/1009_mice_data.xlsx', 1)
mice_data = mice_data[complete.cases(mice_data[, varsx]), ]   # 删除自变量为空的行
mice_data$score_qualified = ifelse(mice_data$score_qualified == 1, 2, 1)  # 置换  score_qualified（1 -> 2, 2 -> 1）
mice_data$id = as.character(mice_data$id)
mice_data[varsy] = lapply(mice_data[varsy], as.numeric)
mice_data[varsx] = lapply(mice_data[varsx], as.factor)
mice_data[varsc] = lapply(mice_data[varsc], as.factor)

# 2.3 gee
fit_t1 = gee_run(mice_data,
        varsx = varsx,
        varsy = varsy,
        varsc = varsc[-5],
        file = "../output/1009_mice_10cv_del_birthweight_c.xlsx")
fit_t2 = gee_run(mice_data,
        varsx = varsx,
        varsy = varsy,
        varsc = varsc[-2],
        file = "../output/1009_mice_10cv_del_week_c.xlsx")
